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Abstract 


Low order equivalent system (LOES) models for the Tu-144 supersonic transport 
aircraft were identified from flight test data. The mathematical models were given in terms 
of transfer functions with a time delay by the military standard MIL-STD-1797A, “Flying 
Qualities of Piloted Aircraft,” and the handling qualities were predicted from the estimated 
transfer function coefficients. The coefficients and the time delay in the transfer functions 
were estimated using a nonlinear equation error formulation in the frequency domain. 

Flight test data from pitch, roll, and yaw frequency sweeps at various flight conditions were 
used for parameter estimation. Flight test results are presented in terms of the estimated 
parameter values, their standard errors, and output fits in the time domain. Data from 
doublet maneuvers at the same flight conditions were used to assess the predictive 
capabilities of the identified models. The identified transfer function models fit the 
measured data well and demonstrated good prediction capabilities. The Tu-144 was 
predicted to be between level 2 and 3 for all longitudinal maneuvers and level 1 for all lateral 
maneuvers. High estimates of the equivalent time delay in the transfer function model 
caused the poor longitudinal rating. 
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1. Introduction 


Flying qualities, or handling qualities as they are also called, are defined as 
“qualities or characteristics that govern the ease and precision with which a pilot is able to 
perform the tasks required in support of an aircraft role” (ref. 1). The handling qualities 
can only be assessed from pilot opinions, but the governing military standards for flying 
qualities offer methods of predicting the handling qualities from estimated transfer function 
coefficients. The Military Specification for Flying Qualities of Piloted Aircraft, MIL-F- 
8785B (ref. 2), was written in the 1960’s for unaugmented aircraft, or aircraft which did not 
have higher-order control systems (HOS). Predicting handling qualities from the open- 
loop transfer function coefficients was an acceptable method for unaugmented, or classical, 
aircraft but the use of complex augmentation systems on aircraft required a different 
mathematical model to describe the aircraft dynamics than the open-loop transfer function. 
Many aircraft built in the 1970’s with high-order control systems were designed without the 
benefit of the specification since the guidelines were not considered to be applicable. A 
revision to the military specification, MIL-F-8785C (ref. 3), was developed in 1980 and was 
the first to recognize augmented aircraft and introduce the low order equivalent systems 
(LOES) concept: 

“The contractor shall define equivalent classical systems which have responses 
most closely matching those of the actual aircraft. Then those numerical 
requirements.. .which are stated in terms of linear system parameters (such as 
frequency, damping ratio and modal phase angles) apply to the parameters of that 
equivalent system rather than to any particular modes of the actual higher-order 
system.” 
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Thus, the military specification suggests that the equivalent system model should have 
parameters which are directly relatable to their classical counterparts. The LOES models are 
linearized, reduced-order models of the actual aircraft response and are used to allow the 
existing flying qualities requirements established for unaugmented, or classical, aircraft to 
be extended to higher order systems. 

The most recent military standard for flying qualities of piloted aircraft, MIL-STD- 
1797A (ref. 4), was written in 1990 and suggests specific LOES models in the assessment 
of flying qualities. The LOES models have the same structure as the classical open-loop 
linear models, but use the pilot input with a pure time delay rather than control surface 
deflections as the model input. 

Mitchell and Hoh (ref. 5) found encouraging results with a number of high-order 
system flight test results and showed in most cases the pitch short-period equivalent 
dynamics are relatable to their unaugmented counterparts. Additionally, a key finding was 
the significant role that the time delay played in the degradation of longitudinal flying 
qualities. Pure time lags produced by the control system directly impacted the pilots’ 
opinions of the handling qualities. 

In terms of the MIL-STD- 1797 A models, a substantial amount of research has been 
performed by Tischler (ref. 6-8) on their identification using frequency response matching. 
This method uses a least squares fit of the Bode plot (magnitude and phase) in the 
frequency domain using the LOES as the model. Though the results of this method have 
been very good, frequency response matching requires substantial data conditioning and 
computation time to estimate accurate parameters. 

A more direct approach was sought by Manning and Gleason (ref. 9) who estimated 
the parameters of the LOES model in the time domain. Time response matching is an 
attractive option since the measured input and output data are all that are required; that is, no 
transformation to the frequency domain is required. 
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The flying qualities of a supersonic transport aircraft may be significantly different 
than current subsonic transports due to the vastly different configurations required for high- 
speed flight. Great concern is placed on the handling of the aircraft during landing 
approach where the aircraft is typically more difficult to control and greater precision is 
required. The NASA High Speed Research (HSR) program, in conjunction with Boeing 
and the Tupolev Design Bureau, performed 19 flight tests of a Russian Tu- 144 supersonic 
transport to establish a new database of information for the development of a U.S. 
supersonic transport in the early 21st century. Data for flying quality and aircraft response 
evaluations were recorded from an array of aircraft configurations and flight conditions. 
From the recorded data, LOES models can be identified, and the flying qualities of the 
aircraft in several different flight regimes can be predicted. The primary objective of this 
research is to predict the flying qualities of the Tu-144 from flight test data. 

This paper begins with a description of the Tu-144 aircraft and the flight test data 
used for flying qualities prediction, as well as a data compatibility check and description of 
all corrections made to the data before analysis. The mathematical models given by the 
military standard are given with all necessary assumptions and the important parameters for 
flying qualities prediction are illustrated. 

Next, the reasoning for the selection of the equation error and output error method 
in the frequency domain for identification of the LOES models is presented. These 
methods are developed completely in reference 10, but a brief development is presented in 
this paper as well. Simulated data was used to validate the mathematical models and the 
identification methods. 

Finally, the methodology was used to identify the model parameters from flight test 
data for twenty-one pitch frequency sweeps, nineteen roll frequency sweeps, and nineteen 
yaw frequency sweeps which were performed at various flight conditions and aircraft 

configurations. The parameters in the identified models were then used to predict the flying 

■ » 

qualities of the Tu-144. The results are presented in terms of the estimated parameters and 
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their standard errors for each maneuver with the predicted handling qualities for each 
maneuver. If the aircraft was not rated level 1, the reason is also given. Suggestions for 
future research are also discussed. The appendices contain a detailed description of the 
flight test maneuver instructions given to the pilots as well as a derivation of the error 
analysis. 



2. Aircraft and Flight Data 


The Tu-144 supersonic transport aircraft is shown in figure 1 and the geometric, 
mass, and inertia characteristics of the aircraft are summarized in table 1. The Tu-144LL 
used for flight tests was a refurbished Tu-144D, tail number 771 14, and re-equipped with 
Kuznetsov NK-321 turbofan engines like those operational on the Russian Air Force 
Tupolev Tu-160 Blackjack bomber. Pitch and roll are controlled with elevons which extend 
along the aft of the entire span of the wing. A rudder on the vertical tail is used for 
directional control. The aircraft is equipped with conventional cable-commanded hydraulic 
actuators which have a parallel electronic input. For stability augmentation, the cable input 
from the pilot is summed with the electronic command from the flight control system. 
Control inputs from the pilots go directly to the actuators via the control cables and do not 
go into the control laws. 

The pitch augmentation system consists of a pitch damper which uses pitch rate 
feedback to provide improved pitch damping of the short period mode. Similarly, there is a 
yaw damper that uses yaw rate feedback. Autopilot and autothrottle were turned off during 
the testing. During takeoff and landing, the nose is drooped for increased visibility and a 
canard above the cockpit is extended for increased stability. The angle of droop for the 
nose for takeoff and landing is 1 1° and 17°, respectively. The canard is extended only for 
stability at lower speeds and is not a control surface. 

For each flight test, the pilot input, control surface deflections, and aircraft responses 
were measured and stored as time histories. Additionally, other variables used to verify the 
flight condition and aircraft configuration were measured. The FDAS variables used for 
data analysis are listed in table 2. Table 3 indicates the coordinate locations of the 
instrumentation measuring aircraft response. 
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All flights of the Tu-144 were performed by Russian pilots at the Zhukovsky Air 
Development Center near Moscow, Russia. Flight tests were performed at speeds ranging 
from Mach 0.3-1 .6 at altitudes of 5,000-50,000 feet and angles of attack ranging from 4°- 
1 1°. The test point number (TPN), test title, aircraft configuration, flight conditions and 
flight number for the maneuvers are listed in table 4. 

For the handling qualities experiment, basic airworthiness sensors on board the Tu- 
144LL were used. No additional instrumentation specific to this experiment was installed. 
The parameters were sampled into analog input channels in a DAMIEN pulse code 
modulation (PCM) data acquisition system. These channels were pre-sample filtered to 
prevent aliasing of data. The filters used were 2-pole low-pass passive RC filters with a 1 
dB per octave rolloff. Attenuation was 3 dB at 200 Hz. Flight test data was taken by 
Tupolev and transferred to the NASA Dryden flight research center where it was made 
available on the flight data access system (FDAS). Further discussion of the 
instrumentation system can be found in reference 11. 

The maneuvers used for handling qualities prediction were frequency sweeps 
performed along the pitch, roll and yaw axes at each flight condition or aircraft 
configuration of interest. A frequency sweep is a commanded oscillation of the controls 
about a trim condition which increases in frequency from the start to finish of the maneuver. 
Frequency sweeps along each axis were chosen since they would excite a wide range of 
frequencies. Additionally, these maneuvers are ideal for frequency response matching since 
they contain a rich spectral content. Doublets were performed along the same axes at each 
flight condition and these data were used to assess the predictive capabilities of the 
estimated models. A doublet is a combination of two pulses of equal amplitude and 
opposite sign in succession. Specifications given to the Russian pilots denoting how the 
frequency sweep and doublet maneuvers were to be flown are taken from reference 12 and 
included in appendix A. 


6 



Translation of the Russian coordinate system to the U.S. coordinate system and the 
calculation of true airspeed from indicated airspeed were both required for the Tu-144. 
These corrections are made at NASA Dryden and the corrected variables are available 
directly from the FDAS system. Corrections to the angle of attack and sideslip angles to 
account for vehicle rotation and corrections to accelerometer measurements due to center of 
gravity offset were not made at NASA Dryden but were completed before data analysis. 

As a part of the data analysis, a data compatibility check was conducted. The 
purpose of this check is to identify and estimate constant offset and scale factor errors in the 
measured response variables due to instrumentation. The error parameters are added to the 
equations of motion and estimated using a maximum likelihood technique (ref. 1 3). If there 
was a consistent and significant effect of these parameters, the data were corrected prior to 
analysis. 

For the longitudinal mode, the outputs were velocity, angle of attack, and Euler pitch 
angle. The lateral outputs were the sideslip angle, Euler roll angle, and Euler yaw angle. 

The main emphasis was on the agreement between the measured angle of attack and 
calculated angle of attack using integrated acceleration measurements for the longitudinal 
mode and the measured and calculated sideslip angle for the lateral mode. Even with the 
linear bias and scale factor error parameters, the data compatibility routine could not 
satisfactorily fit the angle of attack or sideslip angle for a longitudinal or directional 
maneuver, respectively. The resulting fits for angle of attack and sideslip angle, using flight 
data from test point numbers 2.4- 15. IB and 2.4- 15. 3 A, respectively, are shown in figure 2. 

Different parameters, such as the first derivatives of the estimated angles and higher- 
order order terms, were added into the calculations of these two responses to improve the fit, 
but were unsuccessful. The most notable problem with the fit is that the measured data lags 
the estimation from the equations of motion. The cause of this lag is unknown. 
Additionally, the measurement of sideslip angle flattens out at the higher frequencies. This 
may complicate identification of the models for higher frequencies. 
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Another unexplained behavior in the data was the presence of time skews, or 
abnormal time shifts in the data. A time delay exists between the stick input and the control 
surface deflection due to the control system. Once the control surface deflects, however, we 
would expect the aircraft to respond almost immediately. This was not the case in for the 
Tu-144. An example of this type of time skew is shown in figure 3. The elevator deflection 
lags the stick by approximately 0. 1 3 seconds. This lag can be attributed to the control 
system and actuator delays and is normal for a large transport aircraft. The pitch rate then 
lags the elevator deflection by an additional 0.25 seconds. Initially, the time skews were 
believed to be caused by time intervals between sampling different parameters on the 
multiplexed data system; however, this time delay can be at most one time frame long, or 
0.03125 seconds. This accounts for very little of the 0.25 second lag. Another theory is 
that this time delay is real and the aircraft actually responds in this manner. Both Boeing 
and Tupolev have been made aware of the problem, but the cause for the delay is still 
unknown. Similar time skews are apparent with the lateral maneuvers as well. Whatever the 
cause, these time skews are not accounted for in the equations of motion used for data 
compatibility and may be the cause of the poor agreement with the angle of attack and 
sideslip angle estimates. Since there was no information on the magnitude (if any) of the 
time skews, the data was analyzed as recorded and stored on FDAS. 



3. Mathematical Models and Military Standard 

Requirements 


The LOES models given by the military standard are introduced in this section and 
the parameters required for handling qualities prediction are highlighted. In some cases, the 
mathematical models deviate from those suggested by the military standard, but the 
assumptions which lead to the different models are explained. The estimated coefficients 
are directly correlated to the handling qualities criteria, and these criteria for the longitudinal 
and lateral modes for a Class m, Category B or C aircraft are summarized in tables 5-7. In 
every case, the primary goal is to formulate mathematical models which not only describe 
the approximate dynamics of the aircraft, but will also yield the specific parameter estimates 
which lead directly to the handling qualities prediction. 


3.1 Longitudinal Models 

For handling qualities prediction, the military standard places requirements on the 
short period damping ratio (^ sp ), time delay (t 6 ), and the product of the short period natural 
frequency and inverse of the high frequency pitch attitude zero (to sp T 0 2 ). 

The pitch rate and normal acceleration LOES models are given in reference 4 as 
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Note all of the required parameters are found in these models. For a maneuver where the 
velocity is approximately constant, a short period approximation may be made. 

Additionally, when the effect of elevon deflection on the lift is neglected, the normal 
acceleration to stick transfer function may be rearranged to be a function of only parameters 
which appear in the pitch rate to stick transfer function. This is advantageous for parameter 
estimation since it allows the use of another measurement without any additional 
parameters. The new LOES models may be rewritten and expressed with generic 
parameters as 
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B = K e /T 02 , 
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All of the required parameters for handling qualities prediction are present in these 
equations even though the models have been greatly simplified. Equations (3) and (4) 
without the exponential term are identical to the open loop transfer function models with 
elevator deflection, 5 e , as the input. The parameter T e accounts for the time delay between r) e 
and 5 e as well as other possible nonlinearities and added dynamics associated with the 
control system and aircraft augmentation. However, it is important to note that although the 
primary function of x 0 is to account for the time delay between stick and control surface, the 
equivalent system will actually estimate the time delay as the time between stick deflection 
and aircraft response (input to output), or r| e to q, for the longitudinal mode. 

Table 5, given by the military standard, relates the parameter values estimated from 
equations (3) and (4) to the handling qualities of the aircraft. Note the requirement on the 
time delay is 0. 10 for a level 1 aircraft. There has been significant data to suggest that this 
figure is too stringent for large transport aircraft and values of x e of up to 0.4 have still 
resulted in a level 1 pilot rating (ref. 14-16). Nevertheless, the current military standard 
values were used in this report. 

The generic parameters, [k, 1^ A B x e ] T , are introduced to simplify the model and 
ease the workload of the optimizer. Obtaining the military standard transfer function 
coefficients from the generic parameters is a simple algebra problem, but determining the 
errors from the generic parameters required a little more computation. The derivation of the 
error propagation from the generic parameters to the military standard transfer function 
coefficients is given in appendix B. Thus, we can estimate all the necessary parameters for 
handling qualities prediction and their standard errors using equations (3) and (4). 
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3.2 Lateral Models 


For the lateral modes, the military standard requires only 3 parameters to be 
estimated: the roll-mode time constant (T R ), the Dutch roll damping ratio (£ d ), and the 
Dutch roll natural frequency (co d ). The requirements imposed on the values of these 
parameters are indicated in tables 6-7. 

The military standard gives two options For obtaining the Dutch roll damping and 
natural frequency. First, if the ratio of amplitudes of bank angle and sideslip angle 
envelopes in the Dutch roll mode, |<j)/P| d , is large, reference 4 suggests estimating the 
parameters in 


and 


K*s 


s 2 + 2C$©<!)S + (O^ 2 e t p s 


n a 


s + — 

V T sy 


s + — Is 2 + 2^ d (O d s + (O d 2 l 

1r It J 


V 


s + 


T, 


M 


v 


s + 


% 


S -h 




Tla 


V 


s + — 

v T sy 


s + ■ 


l R7 


s 2 + 2^ d (0 d s + co d 2 


(5) 


( 6 ) 


simultaneously. Though this would yield all of the parameters required for handling 
qualities prediction, the models contain a substantial number of parameters which would be 
difficult to estimate accurately. 

If |<(>/P| d is small, the military standard then suggests the use of a second-order 
transfer function relating sideslip angle to rudder pedal deflection: 
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(7) 
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This transfer function is obtained from the Dutch roll approximation and assumes the side 
force due to rudder input is negligible. An alternative transfer function utilizing only the 
Dutch roll approximation relates yaw rate response to yaw control input: 
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Equation (8) adds an additional parameter to the estimation; however, it does offer the ability 
to use the rate measurement instead of the sideslip angle measurement. The estimation of 
the parameters in either (7) or (8) would give the required Dutch roll parameters for 
handling qualities prediction. 

To obtain the maximum roll-mode time constant, reference 4 defines the equivalent 
roll and sideslip transfer functions, respectively, as 
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The military standard does not suggest these two transfer functions be estimated 
simultaneously, but in either (9) or (10), there is still a substantial number of parameters to 
be estimated. If we look only at (9) and assume 


p = <j) or p«s<j). 


then we can rewrite the transfer function as 



K.s 
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( 11 ) 


Additionally, if we assume that the spiral mode time constant, T s , is large, then we can 
simplify (1 1) to 


K<t> 

s 2 + 2^(j,C0(j)S + to,}, 2 
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K^s 2 +Cs + D]e _T P s 
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( 12 ) 


where 

C — 2K ( j,^ ( j ) co ( j ) , 

and 

D^K^ 2 . 

This simplification not only reduces the order of the transfer function, but also eliminates 
one parameter from the estimation. Once again the model is stated in terms of generic 
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parameters as in (3); however, there are still a large number of parameters to be estimated in 
the transfer function. One option to reduce the number of parameters in the estimation is to 
fix the values of those which have already been estimated. For example, if a yaw frequency 
sweep was used to estimate the Dutch roll parameters in equation (7) or (8), then those can 
be fixed in the estimation of ( 12) when using a roll frequency sweep. 

Another alternative which reduces the number of parameters significantly is to 
assume the numerator and denominator quadratics are nearly equal, which assumes the 
aircraft behaves as a first-order system: 


V s + — ' 

Tr 

This transfer function allows only one degree of freedom, and a large modeling error will be 
introduced if this model is used and the assumption of pole-zero cancellation of the 
quadratics is not valid. 
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4. Identification Methodology 


The selected identification techniques were the equation error and output error 
methods in the frequency domain. In the following sections, the reasoning for this choice is 
explained and the selected methods are developed for one and two measurements for the 
longitudinal mode. This development is based on reference 17. The lateral and directional 
modes are easily developed in the same manner. 

4.1 Time and Frequency Domain 

The time domain method used in reference 9 minimized the squared error between 
the measured data and model output for the parameter estimation. Time domain matching 
was an attractive option since the measured output and input data were all that were required. 
However, the optimization technique, when applied to (3) and (4), can have serious 
convergence problems. This can be illustrated by rearranging the transfer function (3) and 
expressing it in the time domain: 

q + k|q + k 0 q = Aq e (t - Tq) + Bq e (t - Te). (14) 

The time delay, t 6 , is a parameter to be estimated; however, any perturbation in x H would 
change the input form. If the optimization algorithm can actually converge to a solution, the 
estimates would likely have high errors. 

On the other hand, rearranging (3) in the frequency domain gives 

-0) 2 q + jcokjq + k 0 q = jcoAfj c e _ -’ (OTe + Bf| e e - j“ Te . 




(15) 
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The most notable benefit to frequency domain analysis is the time delay is an ordinary 
parameter and will not alter the input form during the estimation procedure. Additionally, 
the integration of the state equations is not required for the estimation procedure; the 
mathematics are reduced to algebraic manipulations. The primary disadvantage to frequency 
domain analysis is the data must be transformed for analysis. Frequency domain 
techniques were selected to reduce convergence problems and eliminate numerical 
integration. 

4.2 Equation Error Method 

In the equation error method (EEM), the measured time histories are Fourier 
transformed and the transfer function model is used to match the complex data in the 
frequency domain. Translation to the frequency domain was performed with a high 
accuracy Fourier transform (ref. 1 8) to eliminate translation errors. The typical formulation 
for the equation error method is for only one output measurement. However, the 
formulation can also be extended to two or more output measurements. Both of these 
formulations are briefly presented below, and a full development of the equation error 
method for one output measurement is presented in reference 10. 


4.2.1 One Output Measurement 

In EEM, the sum of squared errors satisfying the equation is minimized. Recall 
equation (15): 


-<*> 2 q + jcok,q + k 0 q = jo)Af| e e J 0510 + Bfi e e - ioyCe . 

When all terms containing unknown parameters are moved to the right-hand side, the 
estimation for the single measurement of pitch rate then becomes: 
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The parameters can then be perturbed until the right-hand side of (16) is equal to the left- 
hand side within an acceptable stopping criterion. Specifically, the cost function to be 
minimized for (16) is 


1 m r 

w=4x 

2 i=i 


.2 

— cof qj +(j©jk, + k 0 )q; -(j©iAe _jft)iTe + Be _j&>lXe jf| e , (17) 


where X represents the vector of estimated parameters, X = [k, Iq, A B x]’ and m represents 
the number of frequencies. A tilde (~) over the variable represents the Fourier transform of 
that variable. Since equation (17) is nonlinear in the parameters, the parameter estimation 
constitutes a nonlinear estimation problem and will require the use of an iterative technique. 
The modified Newton-Raphson technique (ref. 19) was employed for the estimation of the 
parameters. 


4.2.2 Two Output Measurements 

The formulation for two measurements is identical in theory to that for one 

measurement with a few exceptions. If the second measurement of normal acceleration is 
rearranged in the same way as pitch rate, an equation analogous to (16) is formed: 


-to 2 a z = 


-j(oa z - a z 


■— f| e e - ' 0>t e 


(18) 
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The optimization scheme will now have to estimate parameters that satisfy equations (16) 
and (18) simultaneously. The cost function for this method can be written as 



( 19 ) 


where the E subscript denotes experimental, or measured value; that is, 



-0) 2 a, 


( 20 ) 


Then, 


¥(>.) = 




( 21 ) 


is formulated where Y q and Y a ^ are equivalent to the right-hand sides of equations (16) 

and (18), respectively. S vv is the spectral density of the measurement noise estimated from 
the residuals: 
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( 22 ) 


where 


v = Y E -Y (X). 


( 23 ) 
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This weighting matrix is required to account for the different physical values of the pitch 
rate and normal acceleration. 

4.3 Output Error Method 

An alternative approach also developed in reference 10 is the output error method 
(OEM). In this approach, the sum of squared differences between the measured and model 
outputs is minimized. Aside from this difference the development is identical to the 
equation error method. 

4.3.1 One Output Measurement 

The development begins with the transfer function relating pitch rate to stick, 
equation (3): 

q (As + B)e~ TeS 
fi e s 2 + k,s + k 0 

Now the equation is rearranged and converted to the frequency domain leaving the 
measured output, q , segregated on the left-hand side and the remaining terms on the right- 
hand side: 


_ (j( 0 A + B)e- j0)Xe . 

q- 2 , — r-^e- 

—to + jcokj +k 0 


(24) 


Note that the optimization theory is the same; that is, vary the parameter values until the 
right-hand side of (24) becomes acceptably close to the left-hand side. The cost function is 
written as 
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( 25 ) 



where the subscript E once again denotes the experimental value and qj(A,) denotes that 
value of pitch rate obtained from the right hand side of (24). 


4.3.2 Two Output Measurements 

If desired, the measurement of normal acceleration can also be used in OEM. The 
individual equation for normal acceleration analogous to (24) is 


0-Be - j“ te 

a z = 2 • , . He- 

— © + jcok | + kg 


The cost function looks identical to that for EEM: 


(26) 



(27) 


but now Yg is given by 


Y E = [q a z ] T , 


(28) 
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and 



vW=[y, y,J t 


(29) 


is formulated from the right hand sides of equations (24) and (26), respectively. The S vv 
matrix is formed in the same way as in EEM. If desired, the estimates obtained from EEM 
can be used as initial estimates in OEM to decrease convergence time. 


4.4 Output Error Method for Frequency Response Data 

One of the most popular estimation techniques in the frequency domain involves 
frequency response matching (ref. 6-8). This entails a least squares fit of the Bode plot 
(magnitude and phase) in the frequency domain using the transfer function as the model. 
The frequency response from the measured time histories is found from a ratio of the cross- 
spectral density of the input and output to the auto-spectral density of the input. For the 
longitudinal case, this can be written as 


H(co) = 


^t| c q( w ) 

S Tl e Tl e (®) 


(30) 


The accuracy of the model identification depends on the accurate computation of frequency 
response data points from the measured data; subsequently, this requires accurate spectral 
estimates. In order to obtain accurate spectral estimates, reference 7 suggests the use of 
four different data conditioning techniques: digital prefiltering, overlapped/tapered 
windowing, the chirp z-transform, and composite window averaging. These methods not 
only require a significant amount of computation time and effort, but the accuracy of the 
spectral estimates is also a function of the amount of data available. 
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Though frequency response matching is a very common estimation technique in 
parameter estimation, it would be advantageous to use a method which does not require such 
enormous amounts of computational work. 
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5. Analysis of Simulated Data 


Having the model structure predetermined as a LOES by the military standard is 
advantageous since a major step in the system identification process, i.e. model structure 
determination, can be eliminated; however, it can be a great disadvantage if the model is 
inadequate in representing the higher-order system (HOS). As a check on the model 
structure and the estimation algorithms, a simulation case was developed for the longitudinal j 

mode which could judge the performance of the LOES model using outputs generated by a 
higher order system which emulated actual Tu-144 flight dynamics. For the lateral modes, 
simulation cases were developed which checked only the identifiability of the models. Both i 

types of simulations and the conclusions drawn from them are discussed in this section. j 

5.1 Longitudinal Simulation 

The Tu-144 HOS model for the longitudinal mode was created by adding first order 
control system dynamics to the short period mode. The dynamics of the model were 
chosen to approximate the Tu-144 during a maneuver performed at a Mach number of 0.9 1 

at 32,000 feet and angle of attack of b°. Only approximate coefficient values were required 
since the order of the system was of greater importance than the exact model. The HOS 
models were determined to be: 


q _ 12.4s + 8.1 

fie (s + 0.98)(s 2 + 2(0.6)(4.3)s + (4.3) 2 ) 


(31) 
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and 


a*_ = -219.6 

Tie (s + 0.98)(s 2 + 2(0.6)(4.3)s + (4.3) 2 ) 


( 32 ) 


By including the first order control system, an additional pole was added to those of the 
short period mode. 

These models along with a frequency sweep input from measured Tu-144 data were 
used to create simulated time histories of pitch rate and normal acceleration. Gaussian noise 
was added to the simulated outputs to represent random measurement variations. The noise 
level had a standard deviation of 10% of the root mean square value of the simulated output. 
These noisy outputs were then used as “measured” time histories to estimate parameters in 
the LOES models to see how well the LOES estimation of the output could match the 
simulated time histories. Parameters in (3) and (4) were estimated using EEM and OEM 
for both a single output (pitch rate, q) and for two outputs (pitch rate, q, and normal 
acceleration, a z ). The frequency range of interest was 0. 1 rad/sec to 27t rad/sec in 0.0 1 
rad/sec increments. A wider frequency range of 0. 1 rad/sec to 10 rad/sec is suggested by 
the military standard; however, upon spectral analysis of the input signal for all maneuvers 
on the Tu-144, it was found that the signal had little frequency content at frequencies greater 
than 2 tu rad/sec. Thus, the narrower band was selected. Though the estimation is carried 
out in the frequency domain, it is more physically meaningful to compare the estimated 
model in the time domain with the measured output. The time histories of the simulated 
output from (31) and the estimation of the LOES model output from EEM for one output 
are shown in figure 4a. Figure 4b shows the residuals, the difference between simulated 
and model output, for this fit. Even in the presence of measurement and modeling errors, 
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the fit in the time domain is excellent. The application of the other three methods on the 
same data produced similar fits for all outputs. The parameter estimates for the single 
output EEM were 

' K e ] I" 1.972 (0.054)' 

1/T 02 2.048 (0.116) 

£ sp = 0.607(0.018) . 
co sp 2.922(0.061) 

x e J |_0.120(0.004)_ 

The standard errors of the parameters, noted parenthetically next to the estimates, were 6% 

or less. The short period damping ratio, ^ sp , remained approximately the same as in the j 

higher-order model, but the natural frequency, to , decreased by roughly 30%. The 1 

j 

estimates of the static gain, K^, and pitch attitude zero, 1 /Tq 9 , also differed significantly 

from the higher-order system. This indicates some of the effects the addition of the time 
delay and the use of the LOES may have on the parameter estimates. { 

Upon analysis of the normalized pairwise parameter correlation matrices, it was 
discovered that there was a high correlation between k 0 and B for all of the cases except 
EEM for two outputs. The OEM model with one output also had a high correlation 
between k, and A, and the OEM model with two outputs had additional high correlations 
between k, and k,,, A and B and between k,, and A. The correlation matrices for all four 
cases are presented in table 8. A ‘high’ correlation is defined in this paper as one whose j 

absolute value is greater than 0.90. Since the correlation matrix is symmetrical, the upper T 

triangle is blacked out for clarity. 

Insight into the differences between estimation techniques may be evident in the 
correlation matrices. Correlated variables in an estimation are mathematically analogous to 
having more unknowns than equations. There is simply not enough information to find a 
unique solution. Thus, high correlations may reduce the accuracy of the parameter I 
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estimates. In three of the methods, a high correlation occurs between the parameters k 0 and 
B. Together, these two terms represent the static gain of the transfer function, and the model 
is not structured in a way to prevent this correlation. In the single measurement OEM, an 
additional high correlation is present. Finally, for the two measurement OEM, a total of five 
high correlations exist. The data and model structure were held constant for all methods. 
This indicates that the source of the differences between the correlation matrices was the 
estimation technique itself. 

Looking only at the number of high correlations in the estimation, the simulation 
cases illustrated that for the given models and data sets that in this application, the output 
error method was not an appropriate choice for this purpose. Eliminating errors in the 
methodology was critical to the overall quality of the results; thus, the output error method 
was abandoned at this point. EEM for both one and two outputs were still to be used. 

5.2 Lateral Simulation 

The purpose of the lateral simulations was a check on the identifiability of the 
models. A model is said to be identifiable if all of the model parameters can be estimated 
with the identification method and given data. The applicability of a LOES to a high-order 
system was not of primary importance. In these simulations, models given in section 3.2 
were used along with frequency sweep inputs to generate simulated outputs. Noise was 
added as in the longitudinal simulations to the outputs and an attempt was made to estimate 
the parameters in the same model using either EEM or OEM. Attempting to estimate all of 
the parameters at once, using a roll frequency sweep with either (5) and (6) simultaneously 
or (9) individually, or using a yaw frequency sweep with (10), was not successful. The 
parameters were highly dependent on the starting values, especially if the numerator and 
denominator quadratic terms were approximately equal. When these values were not close 
to one another, the model could be estimated, but with substantial computation time and 
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several high pairwise parameter correlations. Another option was to estimate the Dutch roll 
and roll-mode time constant separately using different models. Both (7) and (8) were used 
to estimate the Dutch roll parameters. Each of these models were easily identifiable and the 
parameters all had low standard errors and approximated the known values in the simulation 
within 1%. For the roll-mode time constant, equations (12) and (13) were used. Again, the 
numerator and denominator quadratic terms were important in the estimation, as (12) was 
very difficult to estimate when the quadratics were approximately equal, but the model was 
well-suited when these terms were not nearly equal. Equation (13) was estimated quite 
easily when the simulated data was created with the same model structure, but if a higher- 
order system was used to create the simulated data and the numerator and denominator 
quadratics were not equal, the estimated model had high errors in the fit between simulated 
and estimated output. Thus, (13) is appropriate only if the numerator and denominator 
quadratics are approximately equal. Once again, EEM had fewer high correlations and 
lower standard errors than OEM in all cases. This reiterates the reasons for abandoning 
OEM for this application. 

In short for the lateral modes, the Dutch roll parameters are estimated using either 
(7) or (8) with a yaw frequency sweep. The roll-mode time constant is then estimated using 
a roll frequency sweep and either (12) or (13). If ( 12) is used, the Dutch roll parameters in 
the model are fixed to those estimated from the yaw frequency sweep at the same flight 
condition. 


6. Flight Test Results and Discussion 


This section includes the application of the mathematical models and identification 
methodology to the flight test data. AH maneuvers were analyzed using Matlab v5. 1 on a 
200 MHz Ultra Sparc 2 Unix system. Average time to analyze a single maneuver from raw 
data to final estimate was under two minutes. 


6.1 Longitudinal Results 

Twenty-one longitudinal frequency sweeps were analyzed for the Tu-144. The 
Mach number for these flights ranged from 0.3- 1.6. Parameters of equation (3) and (4) 
were to be estimated using EEM for both one and two outputs (q alone or q and a z , 
respectively) over the same frequency range as the simulations. 

However, EEM for two outputs did not converge for the flight test data. This was an 
indication that there was something inconsistent in the flight test data between the 
measurements of pitch rate and normal acceleration that was not present in the simulations. 
Hence, transfer function coefficients for all 21 sweeps were estimated using the equation 
error method in the frequency domain with one output (q). The estimated transfer function 
coefficients were tabulated along with their standard errors in table 9. Flying qualities were 
then predicted using the estimated transfer function coefficients and MIL-STD- 1 797 A. The 
flying quality predictions are shown in table 10. The longitudinal flying qualities of the Tu- 
144 were predicted to be level 2 or 3 for all maneuvers. The poor ratings were caused by 
high estimates of the time delay in the transfer function. Note that even with the 
modification to relax the time delay mentioned in section 3.1, the Tu-144 would still be a 
level 2 aircraft in all cases. 

Figure 5 contains plots of the values and standard errors of each parameter for every 
test point number. This figure best illustrates the proximity of the results for the repeated 
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maneuvers. For roughly 90% of the estimations, the repeat cases were within 2 a of each 
other. One notable exception is test point number 2.4- 15.1 A. Though the fit of the 
estimation and prediction were acceptable, the values for all of the estimated parameters are 
quite different from the other two repeated maneuvers at the same flight conditions. The 
aircraft configuration and flight conditions were identical in the repeat cases. Since there are 
two cases which contradict the parameter estimates of 2.4-15.1 A, this maneuver would likely 
be discarded. Figures 6-8 summarize the results obtained for a maneuver at Mach 0.3, 0.9, 
and 1.6, respectively. Each of these figures includes a plot of the measured and estimated 
frequency sweep, a plot of the measured and predicted doublet, the flight conditions for the 
sweep and doublet, the parameter values and standard errors, and the flying qualities 
prediction. The predicted doublet was at similar flight conditions and the same aircraft 
configuration. The comparisons are shown in the time domain, although the modeling was 
done in the frequency domain. 

Looking again at figure 5, it can be seen that the estimated parameters for test point 
numbers 2.4-3. 1A,B and 2.4-16. lA,B-2.4-21.1A,B were all comparable. For the first 
parameter value in figure 5, K^, the parameter estimates for similar flight conditions are 
grouped together. These maneuvers were all flown at similar flight conditions and aircraft 
configurations. The extension of the landing gear in test point numbers 2.4-19. 1 A, B had 
little effect on the parameter estimates. The most notable differences in parameter values 
occurred in the high subsonic and supersonic maneuvers, but both the flight conditions and 
aircraft configurations were different in these cases so parameter variations can be expected. 
Neither the Mach number nor the angle of attack varied while the other remained fixed; 
therefore, it cannot be determined how either of these variables individually affect the 
parameter estimates. However, average values can be given over a range of similar flight 
conditions where there was no significant change in the estimated parameters. Table 1 1 lists 
three different ranges of Mach numbers and angles of attack: low subsonic, high subsonic, 
and supersonic with their corresponding range of a and average value of each parameter at 
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that flight condition. Test point number 2.4-15.1 A was discarded for the calculation of the 
averages for the reasons mentioned above. With the exception of K e which had extremely 
small errors, the average values fell within the 2a bounds of nearly every estimate in the 
particular range of flight conditions. The parameter estimates had very low standard errors, 
almost always less than 10%, for all maneuvers. 

In general for all for all of the maneuvers, the estimate of the output for the high 
frequency data fit the measured output better than the estimate at lower frequencies in many 
of the frequency sweep estimations. This was due to the formulation of the equation error 
method which has a weighting on the higher frequencies caused by the to 2 term in the 
equation. The low frequency data was typically in the range of 0.3-0.7 rad/sec, and the 
average estimated short period natural frequency was 1.33 rad/sec, so this was not 
considered a major issue in the accuracy of the estimates. 

However, the inaccuracy of the estimates at lower frequencies did affect the doublet 
predictions. The doublets mostly excited only one low frequency; there was very little 
excitation of the higher frequencies. If the low frequencies are not matched well for the 
frequency sweep, the estimated parameters will not be good predictors for the doublets. 
Conversely, if the sweep matched the low frequencies just as well as the high, then the 
estimated parameters would have good prediction capabilities for any frequency. 

An alternative way of looking at the prediction capabilities is not to use the doublets 
at all, but rather employ the estimated model to predict another frequency sweep at similar 
flight conditions. Figure 9 illustrates the prediction of the frequency sweep from test point 
number 2.4- 16. IB where the model has been identified from test point number 2.4- 16.1 A. 
The sweep is well predicted over a wide range of frequencies. Note that the best prediction 
occurs in the higher frequencies; this reiterates the high frequency weighting mentioned 
above. In general, the larger band of frequencies, the more difficult it will be to estimate a 
model which is good at every frequency. If some a priori knowledge of which frequencies 
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were of interest was given, it would be possible to limit the estimation to a more narrow 
band of frequencies to increase the accuracy in this band. 

The normalized pairwise correlation matrices are not presented for each test 
maneuver, but high correlations appeared with the same parameters as in the simulation: k 0 
and B. Additional high correlations with k, and B were present in the takeoff maneuvers 
(TPN 2.4-3. 1 A, B) and the landing gear extended approach maneuvers (TPN 2.4- 19.1 A, B). 
As mentioned in a previous section, the takeoff configuration differed from the approach 
configuration in the nose position. The flap positions are approximately equal for each of 
these maneuvers. For test point number 2.4-19. 1A,B, the landing gear was extended, and 
even though the parameter estimates were approximately the same as for maneuvers at 
similar flight conditions with the landing gear retracted, the correlations were higher. The 
added high correlations for each of these different maneuvers may be due to the unusual 
configurations at these conditions. 

At all supersonic maneuvers, the standard errors of the parameters were the lowest 
and there were no high correlations. This may not be a function of the airspeed, however, 
but rather the canard, which was extended for all low-speed maneuvers. The canard may 
introduce a nonlinear effect which cannot be accounted for in the linearized model. One test 
of this would be to fly a low speed approach maneuver with the canard retracted and 
investigate whether the errors are present in that estimation. These data were not available. 

6.2 Lateral Results 

For the lateral handling qualities, 19 yaw frequency sweeps and 19 roll frequency 
sweeps were analyzed on the Tu-144. The flight conditions and aircraft configurations were 
identical to those for the longitudinal maneuvers. Once again, EEM for a single output was 
employed. Which output was used depended on the model being identified. For the yaw 
frequency sweeps, equation (7) and (8) were identified over the same frequency range as the 
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simulations; however, the agreement between the measured and estimated time histories of 
the sideslip angle, p, attained from the identification of equation (7) was very poor. This 
was likely due to the poor measurement of P at higher frequencies. The time history of P 
showed no oscillatory motion, and therefore had no frequency content, beyond a rudder 
pedal input frequency greater than roughly n rad/sec. Therefore, the model could not be 
identified at these frequencies and this greatly impacted the parameter estimates. Thus, 
equation (7) was discarded in favor of (8) for the yaw frequency sweeps. Identifying the 
model in equation (8) gave the Dutch roll handling qualities parameters and standard errors 
shown in table 1 2. 

Figure 10 contains plots of the values and 2a error bars for each parameter at every 
test point number to once again show the proximity of the parameter estimates to one 
another. Though many of the repeat cases were within 2 a of one another, it can be seen 
that, in general, there was considerably more scatter in these results. However, the two 
parameters of concern for handling qualities prediction, and tD d , are fairly consistent for 
similar flight conditions. Test point number 2.4- 15. 3 A deviated greatly in relation to the 
other two repeated maneuvers and was therefore discarded. The average parameter values 
for the same ranges of flight conditions used in table 1 1 are shown for the Dutch roll 
parameters in table 13. A greater number of the average values fall outside of the 2a ranges 
of the parameters in the Dutch roll estimation than in the longitudinal estimation. This is 
consistent with the scatter and larger errors associated with the Dutch roll estimation. 
Possible explanations for this scatter are addressed below. Figures 11-13 summarize the 
results obtained for a maneuver obtained at Mach 0.3, 0.9, and 1 .6, respectively. Each of 
these figures contains the same information found in figures 6-8 for the longitudinal 
maneuvers. Note that the agreement between measured and estimated time histories is better 
at higher frequencies. This is once again due to the weighting imposed by EEM to the 
larger values of to. 
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The roll-mode time constant was estimated in two ways. Initially, the first-order roll 
rate to lateral stick transfer function in equation (13) was identified. The estimated 
parameters and standard errors are shown in table 14. Plots of these values with 2 a error 
bars are given in figure 14. Average values of the parameter estimates over the 
aforementioned ranges of flight conditions are listed in table 15. In these estimates, there 
seems to be considerable difference between the parameters at Mach 1.2 and 1.6. It is 
unknown why this difference occurs in these parameters. Other than these, the average 
values all fall fairly close to being within the 2a error bounds of the parameter estimates. 

Figures 15-17 summarize the results obtained for a maneuver obtained at Mach 0.3, 
0.9, and 1 .6, respectively. Each of these figures contains the same information found in 
figures 6-8 for the longitudinal maneuvers. Note that for the predicted time histories, the 
model fits the data as well as can be expected from a first-order model, but misses some of 
the higher-order dynamics in the measured data. These dynamics may be caused by either 
the rigid body motion of the aircraft or twisting and bending of the wing and body. This 
gives some indication that the assumption of pole-zero cancellation of the quadratic terms 
mentioned in section 3.2 may not be valid. The agreement between the measured and 
estimated time histories are quite good, but do have a slight mismatch in amplitude at higher 
frequencies (to > Jt rad/sec). The measured and estimated data are in phase with one 
another, so the mismatch can likely be attributed to a problem in the estimation of the static 
gain. Note that the formulation of equation (13) in EEM would not have the same CO 2 
weighting on the frequencies that was present in previous two estimations. 

The roll-mode time constant was also estimated using a hybrid of equation (12). 

The term ‘hybrid’ is used because, although using the form of equation (12) is used, the 
values of the Dutch roll parameters, and (0 d , are fixed during the estimation to values 
obtained from the yaw frequency sweep analysis. This allows a higher-order model to be 
used for the identification, but eliminates two of the parameters to be estimated. The 
estimated parameters and standard errors are shown in table 16. Plots of these values with 
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2a error bars are given in figure 18. Average values of the parameter estimates over the 
aforementioned ranges of flight conditions are listed in table 17. There was some scatter in 
the parameter estimates once again, especially between the Mach 1.2 and 1.6 maneuvers 
again, but the most interesting results were in the estimation of the roll-mode time constant. 
For all three flight condition ranges, the roll-mode time constant was fairly consistent with 
the first-order estimation. For a Mach number of 0.3 to 0.4 and angle of attack of 8.2° to 
1 1°, the value of the T R varied by only 0.017 seconds. The differences for the second and 
third flight condition ranges were 0.04 seconds and 0.12 seconds, respectively. Figures 19- 
21 summarize the results obtained for a maneuver obtained at Mach 0.3, 0.9, and 1 .6, 
respectively. Each of these figures show the fit, parameter estimates, and handling qualities 
prediction for a lateral/directional maneuver in a format similar to that of figures 6-8 for the 
longitudinal maneuvers. Note that for the predicted time histories, the modeftries to fit 
some of the higher-order dynamics in the measured data and fits the measured data slightly 
better than the first-order model. However, the agreement between the measured and 
estimated time histories is not as good as the first-order model in terms of the root-mean- 
square of the difference between the measured and estimated output. In this estimation, 
similar mismatches occurred at higher frequencies (0) > K rad/sec), but now the estimated 
amplitude was much higher than the measured. 

Note once again that equation (12) was formulated in EEM without a weighting on 
the higher frequencies. The Dutch roll parameters used as fixed values in the model, 
however, were estimated from a model which did weight the higher frequencies. This 
dichotomy may be the cause for the poor agreement at higher frequencies in these estimates 
Nevertheless, the phase at higher frequencies was still identical for the measured and 
estimated time histories, and the amplitude mismatch at high frequencies would not likely 
affect the estimate of the roll-mode time constant. 

With the Dutch roll and roll-mode time constant parameters estimated, a flying 
qualities prediction was made. These parameters are tabulated along with the predicted 
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flying qualities in table 18. Note that the values listed in the table for the roll-mode time 
constant are from the first-order model estimates, but using the third-order hybrid model 
estimates would not affect the handling qualities prediction. The Tu-144 was predicted to 
be level 1 for all lateral maneuvers. Another important note, however, is that the military 
standard places no requirements on the time delay for the lateral handling qualities. 

Additional sources of error are present in the lateral estimations that were not seen in 
the longitudinal mode. The largest probable source of error and the likely cause of the 
scatter in the estimates is due to the coupling of the roll and yaw motion in the aircraft. This 
is attributed not only to the aerodynamic cross-derivatives, but also to the control system. 

For example, when the pilot deflects the rudder pedal, the aileron moves along with the 
rudder. This introduces roll dynamics to the aircraft which may not be accounted for in the 
assumptions made for the yaw rate to rudder pedal transfer function. Though a source of 
error, this coupling will not likely affect the handling qualities prediction. 

The pairwise correlation matrices for the yaw frequency sweep analysis showed no 
high correlations between any of the parameters. For the first-order roll frequency sweep 
transfer functions, there was a high correlation between K p and 1/T R . The third-order 
hybrid transfer function had high correlations between K^, and D. Recall that the mismatch 
at the higher frequencies for the first-order roll model output time histories was disregarded. 
This was because the cause was likely a problem with the static gain; however, it is now seen 
that the static gain is correlated with 1/T R , the parameter of interest. Though normally this 
would be of concern, the third-order hybrid model did not have a high correlation associated 
with the estimation of 1/T R ; the parameters and D are both numerator terms. 

Additionally, the time constants were approximately equal for either the first or third-order 
hybrid model; thus, performing both estimations gives more confidence in the handling 
qualities prediction. This also indicates that even if the assumption of pole-zero cancellation 

in the quadratic terms of the roll transfer function was incorrect, it had little effect on the 

■* 

estimate of T r . -- •• 
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7. Conclusions and Suggestions for Future Work 


The equation error method was applied in the frequency domain for the estimation 
of parameters in the specified transfer functions. These transfer functions were low-order 
equivalent system (LOES) models of the aircraft for pitching, rolling, and yawing motion. 
Some of the transfer functions were altered from those specified in the governing military 
standard to create models that could be more accurately estimated. The change in the 
models did not sacrifice the physical meaning of the estimated handling qualities 
parameters. Simulation cases were developed which emulated the Tu-144 supersonic 
transport and the LOES models were identified using the specified estimation technique. 
The agreement of the simulated data with the identified model response was very good. The 
same estimation procedure was applied to 21 pitch frequency sweeps, 19 yaw frequency 
sweeps, and 19 roll frequency sweeps from flight tests of the Tu-144. Parameter estimates 
and their standard errors were typically comparable between maneuvers repeated at the same 
flight conditions, and the agreement between measured and estimated output in the time 
domain was excellent for all maneuvers. The Tu-144 was predicted to have level 2 or 3 
handling qualities in the longitudinal mode and level 1 in the lateral modes. The reason for 
the poor longitudinal rating was high estimates of the time delay in the transfer function. 

In general, the new formulation of the military standard transfer functions and the 
use of the equation error method in the frequency domain provided an excellent method of 
estimating the handling qualities parameters. The parameter values were estimated with 
standard errors typically less than 10% and could be estimated from the raw data in less 
than two minutes of computer time. Additionally, the first-order approximation of the roll 
rate to lateral stick transfer function yielded nearly identical estimates of the roll-mode time 
constant as the higher-order model which required a priori knowledge of the Dutch roll 
parameters. 
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In the future, the most obvious next step would be to compare the predictions to 
actual pilot ratings of the Tu-144. One major question would be how the pilots rate the 
aircraft in longitudinal motion, and if the results support the notion that the military standard 
requirement for time delay on large transport aircraft is too stringent. 

There are no requirements on the types of maneuvers that must be flown, but only 
that the maneuvers properly excite the modes to be analyzed. Different maneuver forms 
(such as 3-2- 1-1 or 2-1-1 pulses) can be designed which minimize flight time required to 
perform the maneuver while still offering enough information in the data. Another 
maneuver would be to perform step inputs along the roll axis, estimate the roll-mode time 
constant directly from the time history, and compare the results to those obtained from the 
frequency sweep. Maneuvers could also be performed which segregated the parameter 
variations with either Mach number or angle of attack. 

Additionally, frequency response matching is one of the more popular methods of 
parameter estimation when the input form is a frequency sweep. A comparison of the 
results obtained in this paper and those obtained from frequency response matching would 
contribute to the field of parameter identification. 
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Appendix A - Maneuver Description 


The pilots and test engineers were briefed before the flight test program and the 
experiment objectives and test procedures were explained. Descriptions of the individual 
test procedures for each of the maneuvers used in the handling qualities prediction are taken 
from reference 10 and given below: 

Pitch Frequency Sweep 

Conditions: Minimal Turbulence, Autopilot Off, Autothrottle Off, No Throttle Position 

Changes 

1 . A test engineer in the back of the airplane should carefully monitor the control 
column position, load factor, and pitch attitude time history responses plotted on a 
computer screen in real time during this maneuver. The time scale should be 
expanded so that a frequency limit of TBD cycles per second for the input can be 
easily judged. Once this frequency limit is reached, the engineer should tell the pilot 
to stop making inputs to prevent excitation of the 2.5 cycle per second first bending 
mode of the airplane. 

2. Trim airplane for hands off level flight. Do not retrim during the remainder of this 
maneuver. 

3. Begin data recording and record 5 seconds of hands off level flight data. 

4. Slowly cycle the control column back and forth with an amplitude large enough to 
obtain +/- 0.2g load factor and/or +/- 5 to 15 degree pitch attitude excursions. 

Make two complete 20 second cycles of the control for a total of 40 seconds of 
input. The cycling of the control should be centered around a position that produces 
airplane oscillations that center around the trim pitch attitude. Control wheel and 
rudder pedals should be used to minimize roll and yaw response. 
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5. Slowly increase the frequency of the input. Adjust input amplitude so that the 
amplitude of airplane motion remains about the same as in step 4. It is important to 
increase frequency slowly, so that there is enough middle frequency content in the 
data. When the amplitude of airplane motion drops off sharply or the input 
frequency reaches TBD cycles per second, stop the input. Should a structural mode 
become excited, terminate input immediately. The combined inputs for steps 4. and 
5. should last about 80 - 100 seconds. 

6. Record 5 seconds of hands off data at the end of the maneuver. 

Roll Frequency Sweep 

Conditions: Minimal Turbulence, Autopilot Off, Autothrottle Off, No Throttle Position 

Changes 

1 . A test engineer in the back of the airplane should carefully monitor the control wheel 
position, roll rate, and roll angle time history responses of the airplane plotted on a 
computer screen in real time during this test. The time scale should be expanded so 
that a frequency limit of TBD cycles per second for the input can be easily judged. 
Once this frequency limit is reached, the engineer should tell the pilot to stop 
making inputs to prevent excitation of the 2.5 cycle per second first bending mode 
of the airplane. 

2. Trim airplane for hands off level flight. Do not retrim during the remainder of this 
maneuver. 

3. Begin data recording and record 5 seconds of hands off level flight data. 

4. Slowly cycle the control wheel back and forth with amplitude large enough to obtain 
+/- 5 to 15 degree roll angle excursions. Make two complete 20 second cycles of 
the control for a total of 40 seconds of input. The cycling of the control should be 
centered around a position that produces airplane oscillations that center around a 
wings level roll angle. Rudder pedals should only be used if the airplane 
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oscillations do not remain centered about the initial heading angle. Control column 
should be used to minimize pitch response. 

5. Slowly increase the frequency of the input. Adjust input amplitude so that the 
amplitude of airplane motion remains about the same as in step 4. It is important to 
increase frequency slowly, so that there is enough middle frequency content in the 
data. When the amplitude of airplane motion drops off sharply or the input 
frequency reaches TBD cycles per second, stop the input. Should a structural mode 
become excited, terminate input immediately. The combined inputs for steps 4. and 
5. should last about 80-100 seconds. 

6. Record 5 seconds of hands off data at the end of the maneuver. 

Yaw Frequency Sw eep 

Conditions: Minimal Turbulence, Autopilot Off, Autothrottle Off, No Throttle Position 

Changes 

1 . A test engineer in the back of the airplane should carefully monitor the rudder pedal 
position, yaw rate, and heading angle time history responses of the airplane plotted 
on a computer screen in real time during this test. The time scale should be 
expanded so that a frequency limit of TBD cycles per second for the input can be 
easily judged. Once this frequency limit is reached, the engineer should tell the pilot 
to stop making inputs to prevent excitation the 2.5 cycle per second first bending 
mode of the airplane. 

2. Trim airplane for hands off level flight. Do not retrim during the remainder of this 
maneuver. 

3. Begin data recording and record 5 seconds of hands off level flight data. 

4. Slowly cycle the rudder pedals back and forth with an amplitude large enough to 
obtain +/- 5 to 15 degree heading angle excursions. Make two complete 20 second 
cycles of the control for a total of 40 seconds of input. The cycling of the control 
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should be centered around a position that produces airplane oscillations that center 
around the initial heading angle. Control wheel should only be used if the airplane 
oscillations do not remain centered about a wings level roll angle. Control column 
should be used to minimize pitch response. 

5. Slowly increase the frequency of the input. Adjust input amplitude so that the 
amplitude of airplane motion remains about the same as in step 4. It is important to 
increase frequency slowly, so that there is enough middle frequency content in the 
data. When the amplitude of airplane motion drops off sharply or the input 
frequency reaches TBD cycles per second, stop the input. Should a structural mode 
become excited, terminate input immediately. The combined inputs for steps 4. and 
5. should last about 80-100 seconds. 

6. Record 5 seconds of hands off data at the end of the maneuver. 

Pitch Doublet 

Conditions: Minimal Turbulence, Autopilot Off, Autothrottle Off, No Throttle Position 

Changes 

1 . Trim airplane for hands off level flight. Do not retrim during the remainder of this 
maneuver. 

2. Begin data recording and record 5 seconds of hands off level flight data. 

3. Pull back on control column sharply and hold input for 5 seconds, push forward on 
control column sharply and hold input for 5 seconds, and then release the control 
column to neutral position. Inputs should be large enough to produce +/- 0.2g load 
factor and/or +/- 5 to 1 5 degree pitch attitude excursions. Control wheel and rudder 
pedals should be used to minimize roll and yaw response. 

4. Record 60 seconds of hands off data at the end of the maneuver. 


44 



Roll Doublet 

Conditions: Minimal Turbulence, Autopilot Off, Autothrottle Off, No Throttle Position 

Changes 

1 . Trim airplane for hands off level flight. Do not retrim during the remainder of this 
maneuver. 

2. Begin data recording and record 5 seconds of hands off level flight data. 

3. Rotate control wheel sharply one direction and hold input for 5 seconds, rotate 
control wheel sharply the other direction and hold input for 5 seconds, and then 
release the control wheel to neutral position. Inputs should be large enough to 
produce +/- 5 to 15 degree roll angle excursions. Control column should be used to 
minimize pitch response. 

4. Record 60 seconds of hands off data at the end of the maneuver. 

Yaw Doublet 

Conditions: Minimal Turbulence, Autopilot Off, Autothrottle Off, No Throttle Position 

Changes 

1 . Trim airplane for hands off level flight. Do not retrim during the remainder of this 
maneuver. 

2. Begin data recording and record 5 seconds of hands off level flight data. 

3. Push one rudder pedal sharply and hold input for 5 seconds, push the other rudder 
pedal sharply and hold input for 5 seconds, and then release the rudder pedals to 
neutral position. Inputs should be large enough to produce +/- 5 to 15 degree 
heading angle excursions. Control column should be used to minimize pitch 
response. 

4. Record 60 seconds of hands off data at the end of the maneuver. 
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Appendix B - Error Propagation 


The standard errors of the generic parameters introduced in (3) and (4) are not 
equivalent to those of the transfer function coefficients of the models from the military 
standard; however, the estimated standard errors can be used to estimate the errors of the 
desired coefficients through a linearized error propagation formula. Reference 20 explains 
uncertainty analysis in detail, but a short development of the theory and application to an 
example on the longitudinal mode is presented below. 

Consider a general case in which an experimental result, x, is a function of N 
variables, y,: 

x = x(y,,y 2 , ... ,y N ). (33) 

Equation (33) defines how to determine x from the known value of the variables y ( . The 
uncertainty in the result is given as 
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where the U xi represent the uncertainties in the dependent variables y r 

The uncertainties of the estimated parameters and their functional dependence on the 
military standard coefficients are known; thus, equation (34) can be used to determine the 
desired standard errors. For example, the short period damping ratio is given in terms of 
the estimated generic parameters as 
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The errors of the final two parameters, K H and x, are estimated directly. 

These error propagation formulae were used for all of the estimated parameters; 
however, they are only linear approximations of the errors. The first-order analytical 
functions were validated using a Monte Carlo simulation which multiplied a random variable 
with a mean of zero and variance of one by the standard error of each estimated parameter. 
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The product was then scaled by adding the estimated parameter value. For example, for the 
parameter a new parameter, k^’, was created as 

k 0 = k 0 +a(k 0 )*r, (39) 

where r represents a Gaussian distributed random variable. The Monte Carlo estimate of 
the short period natural frequency is then 

( 4 °) 

A random number generator was used to create 1000 values for r and thus 1000 Gaussian 
distributed values for co sp ’. The mean and standard deviation of the Monte Carlo simulation 
should be very close to that estimated by the analytical functions. The other variables were 
checked in the same fashion. Figure 22 illustrates the distribution of each of the checked 
variables for the parameters estimated from test point number 2.4-16.1 A. A small table on 
each plot indicates the mean value and standard deviation for both the analytic function and 
the Monte Carlo simulation. These values indicate that the simulation validated the 
analytical functions. 


48 


iiiiliiiiiililiiiiiiiiiiin ii 1 1 lit' inn n hi mi si 


Table 1 : Summary of geometric, mass, and inertia characteristics of the Tu- 144. 


Length 

196 ft 10 in (60.0 m) 

Span 

88 ft 7 in (27.0 m) 

Nose Tip to Leading Edge of MAC 

98 ft 8 in (30. 1 m) 

Length of MAC 

76 ft 5 in (23.3 m) 

Wing Area 

4716 ft J (438 m 2 ) 

Wing Aspect Ratio 

1.66 

Wing Sweep, Inboard Portions 

76 deg 

Wing Sweep, Main Panels 

57 deg 

Weight* 

303,000 lb (138,000 kg) 

Roll Axis Moment of Inertia, I xx * 

38,805,000 lbf-ft 2 (1,635,000 kg-m 2 ) 

Pitch Axis Moment of Inertia, I yy * 

4 1 7,797,000 lbf-ft 2 ( 1 7,606,000 kg-nL) 

Yaw Axis Moment of Inertia, I Z7 * 

450,222,000 lbf-ft' (18,973,000 kg-nr) 

Roll-Yaw Product of Inertia, I xz * 

-6,486,000 lbf-ft' (-273,000 kg-m") 


* Average values over all yaw frequency sweeps. 
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Table 2: Measured parameters used for data analysis. 


Rate 1 

32 Hz 

32 Hz 

64 Hz 

64 Hz 

N 

X 

<N 

cr, 

64 Hz 

32 Hz 

32 Hz 

16 Hz 

32 Hz 

N5 

X 

<N 

PT. 

32 Hz 

32 Hz 

64 Hz 

64 Hz 

16 Hz 

16 Hz 

16 Hz 

32 Hz 

N 

X 

(N 

rr t 

i Sparse 

8 Hz 

Sparse 

Sparse 

Sparse 

Sparse 

Sparse 

Sparse f 

Accuracy 

±1.5% 

m 

+i 

$ 

*n 

Ti 

±1.5% 

±1.5% 

±1.5% 

$ 

+1 

±1.2% 

±1.2% 

±2.0% 

q 

ri 

■H 

±2.0% 

±2.5% 

±2.5% 

±2.5% 


vq 

o 

+i 


1±0.5% 

±0.5% 


±4.0% 







Range 

-100- +250 mm 

si 

o 

GO 

+1 

±125 mm 

SI 

■o 

+1 

Si. 

V 

"O 

in 

fN 

+1 

SI] 

u 

"O 

in 

CN 

+i 

±90 deg 

Si 

<L> 

■o 

o 

O' 

+1 

±180 deg 

±18 deg/sec 

±6 deg/sec 

CJ 

t/a 

sxj 

o> 

vO 

+1 

±0.5 g’s 

±0.5 g’s 

C/3 

'si 

q 

cn 

+ 

» 

q 


0 - 22000 m 

0.0 - 2.5 

±25 deg 

si 

"O 

m, 

Ti 


0-100 metric tons 







+ Sign Convention 

Pull 

Clockwise 

Right Pedal 

Trailing Edge Down 

Trailing Edge Down 

Trailing Edge Left 

c 

£ 

o 

Q 

Si. 

c 

£ 

-C 

si 

5 

Nose Up 

Clockwise from North 

Right Wing Down 

Nose Up 

Clockwise 

Forward 

S 

sc 

5 

a 

D 

Always + 

a 

D 

Always + 

Nose Up Relative to Flight Path 

Nose Left Relative to Flight Path 

Always + 

Always + 

Always + 

a 

O 

Always + 

Always + 

Always + 

i 

Units 

UI HI 

Si 

■8 

111 m 

Gl] 

-8 



■S 

Si. 

-8 

Dl 

•8 

deg/sec 

§ 

~5l 

8 

c/a 

U 

T3 

C/3 

"si 

05 

r Sl 

05 

'si 

knots 

feet 

No Units 


si 

■8 


metric tons 

% MAC 

E 

"s 

si 

r E 

Ol 

"e 

si 

E 

si 

Parameter Description 

Control Column Position 

Control Wheel Position 

Rudder Pedal Position 

Left Elevon 1 Position 

Right Elevon 1 Position 

Upper Rudder Position 

Roll Euler Angle 

Pitch Euler Angle 

Heading Angle 

Body Axis Roll Rate 

Body Axis Pitch Rate 

Body Axis Yaw Rate 

Longitudinal Acceleration 

Lateral Acceleration 

Vertical Acceleration 

True Airspeed 

Geopotential Pressure Altitude 

Mach Number 

Angle of Attack 

Sideslip Angle 

Post Flight Computed Gross 
Weight 

Total Fuel Quantity Remaining 

Post Flight Computed Long. CG 

Post Flight Computed Vertical CG 

Roll Axis Moment of Inertia 

Pitch Axis Moment of Inertia 

Yaw Axis Moment of Inertia 

Roll - Yaw Product of Inertia 
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1 PPED 

| ERL 
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1 RUDDER 
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I IZZ 

I 1X7. 
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Table 3: Coordinate location of instrumentation. 


Parameters 

X (ft) 

Y (ft) 

Z (ft) 

Airspeeds & Pressures 

-4.4 

-2.7 

0 

Angle of Attack 

+ 15.7 

-1.4 

-2.4 

Sideslip 

+ 15.7 

+ 1.0 

0 

Rate Gyros & Accelerometers 

+ 106.6 

-1.3 

+2.46 


Notes: 

• The origin of the coordinate system is at the base of the nose boom/tip of nose cone. 

• Axis System: 

+X is measured longitudinally from nose to tail 
+Y is measured vertically up 
+Z is measured laterally out the right wingtip 
Note that this is a left-handed coordinate system. 
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Flight Conditions 
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Takeoff Pitch Frequency Sweep 
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Takeoff Roll Frequency Sweep i 
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52 


i ii i. ii i ii ii mill 




Table 4: Continued. 
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Table 4: Continued. 
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Table 4: Concluded. 
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Table 5: Recommended Short Period Requirements for Class III, 
(a) Category B and (b) Category C 


(a) Category B 


Level 

Min ; sp 

■BBS 1 

x, sec 

Min ® SP T 02 

1 

■kesbh 

2.00 

0.10 

1.00 1 

2 


2.00 


0.60 

3 




- 


(b) Category C 


Level 

Min ^ sp 

Max ; sp 

x, sec 

Min C0 sp T e2 

1 


1.30 

■ . 

1.40 

2 


2.00 

■ 

0.70 

3 

- 

- 


- 


Table 6: Recommended Roll-Mode Time Constant for Class III, Category B and C 


Level 

Max T r , sec 

1 

1.4 

2 

3.0 

3 

10.0 


Table 7: Recommended Dutch Roll Frequency and Damping for Class III, 
(a) Category B and (b) Category C 


(a) Category B 


Level 

Min £ d 

Min C/o d 

Min © d 

I 

0.08 


0.4 

2 

0.02 

0.05 

0.4 

3 

0 

— 

0.4 


(b) Category C 


Level 

1 KSSKS 

Min £ d (0 d 

Min co d 

1 

0.08 

0.10 

0.4 

2 

0.02 

0.05 

0.4 

3 

0 

— 

0.4 i 
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Table 8: Normalized pairwise parameter correlation matrices from simulated data for (a) 
EEM for one measurement, (b) EEM for two measurements, (c) OEM for one 
measurement, and (d) OEM for two measurements. Shaded values indicate those defined as 

having a high correlation. 



(a) 



(b) 



(c) 



(d) 
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Table 9: Summary of longitudinal parameter estimates. 
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* Parentheses Denote Standard Errors 



Table 9: Concluded. 
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|: Parentheses Denote Standard Errors 












































































































Table 10: Summary of longitudinal handling qualities predictions. 
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.4-16.1 B Canard Ext. Effect Pitch 0 774 0.200 1.612 B 3 Te > 0.20 

Frequency Sweep . 

Parentheses Denote Standard Errors 













































Table 10: Concluded. 
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Parentheses Denote Standard Errors 



Table 1 1 : Average parameter estimates for different flight conditions. Pitch rate to 
longitudinal stick transfer function coefficients. 


Flight Condition 

K e 

1 / T e 2 

Csp 

^sp 

^0 

0.3 < M < 0.4 
8.2 <a< 11.0 

1.2 

0.74 

0.76 

1.15 

0.20 

0.8 < M < 0.9 
6.0 < a < 6.2 

1.9 

1.2 

0.64 

2.2 

0.22 

1.2<M< 1.6 
4.8 < a < 6.4 

0.75 

0.39 

0.36 

1.55 

0.19 
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Table 12: Summary of Dutch roll parameter estimates. 
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Table 13: Average parameter estimates for different flight conditions. Yaw rate to rudder 

pedal input transfer function coefficients. 


Flight Condition 
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0) d 


0.3 < M < 0.4 
8.2 < a < 1 1.0 

0.36 

0.47 

0.31 

0.83 

0.15 

0.8 < M < 0.9 
6.0 < a < 6.2 
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1.2 < M < 1.6 
4.8 < a < 6.4 
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Table 15: Average parameter estimates for different flight conditions. Roll rate to lateral 
stick input transfer function coefficients. First-order model. 


Flight Condition 

K P 

1/T r 

X P 

0.3 < M < 0.4 
8.2 < a< 1 1.0 

8.4 

2.4 

0.10 j 

0.8 < M < 0.9 
6.0 < a < 6.2 

16.6 

3.2 

0.09 j 

1.2 < M < 1.6 
4.8 < a < 6.4 

6.7 

2.0 

0.14 



















Table 16: Summary of roll-mode time constant estimates. Third-order hybrid model. 
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* Parentheses Denote Standard Errors 














































































































Table 16: Concluded. 
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Table 17: Average parameter estimates for different flight conditions. Roll rate to lateral 
stick input transfer function coefficients. Third-order hybrid model. 


Flight Condition 

-©- 


®<j> 

I/Tr 

T P 

0.3 < M < 0.4 
8.2 < a < 1 1.0 

11.6 

0.18 

0.72 

2.5 

0.12 

0.8 < M < 0.9 
6.0 < a < 6.2 

23.8 

0.30 

0.83 

2.8 

0.13 

1.2 <M< 1.6 
4.8 < a < 6.4 

8.7 

0.16- 

1.5 

2.6 

0.17 


























Table 18: Summary of lateral handling qualities predictions. Roll-mode time constant 

from first-order model shown. 
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* Parentheses Denote Standard Errors 

Note that “2/3” in TPN indicates the use of both a yaw and a roll frequency sweep. 
























































Table 18: Concluded. 
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Figure 1: Three-view drawing of the Tu-1 44 aircraft. 
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(a) 



Figure 2: Data compatibility results comparing measured and calculated (a) angle of attack 

and (b) sideslip angle. 
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Figure 3: Time skew in response of aircraft to control surface deflection. Longitudinal 

doublet, test point number 2.4-15.4A. 
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Figure 4: (a) Simulated and model time histories for pitch rate for EEM, and (b) residuals 

from output Fit, 
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Figure 5: Summary of parameter estimates with 2 o error bars for pitch rate to longitudinal 

stick transfer function. 








Figure 6: Comparison of measured, estimated, and predicted time histories of pitch rate. 
Test point number (a) 2.4- 1 6. 1 A for estimation, and (b) 2.4- 1 6.4 A for prediction. 
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Figure 7: Comparison of measured, estimated, and predicted time histories of pitch rate. 
Test point number (a) 2.4-15. IB for estimation, and (b) 2.4-15.4B for prediction. 
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Figure 8: Comparison of measured, estimated, and predicted time histories of pitch rate. 
Test point number (a) 2.4-12.1B for estimation, and (b) 2.4-12.4B for prediction. 
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Figure 9: Comparison of measured and predicted time histories of pitch rate. Test point 
number 2.4- 16. IB predicted with 2.4- 16.1 A parameter estimates. 
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Figure 10: Summary of parameter estimates with 2a error bars for yaw rate to rudder pedal 

transfer function. 
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Figure 10: Concluded. 
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Figure 11: Comparison of measured, estimated, and predicted time histories of yaw rate. 
Test point number (a) 2.4-1 6. 3A for estimation, and (b) 2.4-16.6A for prediction. 
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Figure 12: Comparison of measured, estimated, and predicted time histories of yaw rate. 
Test point number (a) 2.4- 15. 3 A for estimation, and (b) 2.4-1 5. 6B for prediction. 
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Figure 13: Comparison of measured, estimated, and predicted time histories of yaw rate. 
Test point number (a) 2.4-12.3A for estimation, and (b) 2.4-12.6A for prediction. 
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Figure 14: Summary of parameter estimates with 2cr error bars for roll rate to lateral stick 

transfer function. First-order model. 
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Figure 14: Concluded 
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Figure 15: Comparison of measured, estimated, and predicted time histories of roll rate. 
First-order model. Test point number (a) 2.4-1 7. 2A for estimation, and (b) 2.4-17.5A for 

prediction. 
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Figure 16: Comparison of measured, estimated, and predicted time histories of roll rate. 
First-order model. Test point number (a) 2.4-15. 2A for estimation, and (b) 2.4- 15.5B for 

prediction. 
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Figure 17: Comparison of measured, estimated, and predicted time histories of roll rate. 
First-order model. Test point number (a) 2.4-12.2B for estimation, and (b) 2.4-12.5B for 

prediction. 
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Figure 1 8: Summary of parameter estimates with 2 a error bars for roll rate to lateral stick 
transfer function. Third-order hybrid model. 
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Figure 18: Concluded 
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Figure 19: Comparison of measured, estimated, and predicted time histories of roll rate. 
Third-order hybrid model. Test point number (a) 2.4-17.2A for estimation, and (b) 2.4- 

17. 5 A for prediction. 
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Figure 20: Comparison of measured, estimated, and predicted time histories of roll rate. 
Third-order hybrid model. Test point number (a) 2.4- 15.2A for estimation, and (b) 2.4- 

15. 5B for prediction. 
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Figure 2 1 : Comparison of measured, estimated, and predicted time histories of roll rate. 
Third-order hybrid model. Test point number (a) 2.4-1 2. 2B for estimation, and (b) 2.4- 

12.5B for prediction. 
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